library(magrittr)
library(tidyverse)
library(Seurat)
library(readxl)
library(cowplot)
library(colorblindr)
library(viridis)
library(magick, lib.loc = "/home/uhlitzf/miniconda3/lib/R/library")
library(ggpubr)

theme_cowplot2 <- function(...) {
  theme_cowplot(font_size = 16, font_family = "sans", ...) %+replace%
    theme(strip.background = element_blank(),
          panel.background = element_rect(fill = "transparent", color = NA),
          plot.background = element_rect(fill = "transparent", color = NA),
          panel.border = element_blank())
}
theme_set(theme_cowplot2())

remove_xaxis <- theme(axis.title.x = element_blank(),
                      axis.text.x = element_blank(),
                      axis.ticks.x = element_blank(),
                      axis.line.x = element_blank())

remove_guides <- guides(color = F, fill = F, shape = F, alpha = F)
### load all data ---------------------------------

helper_f <- function(x) ifelse(is.na(x), "", x)
markers_v5 <- yaml::read_yaml("/home/uhlitzf/spectrum_tme/_data/small/signatures/hgsc_v6_major.yaml")

helper_f2 <- function(x) select(unnest(enframe(x, "subtype", "gene"), cols = gene), gene, subtype)
markers_v5_super <- lapply(yaml::read_yaml("/home/uhlitzf/spectrum_tme/_data/small/signatures/hgsc_v6_super.yaml"), helper_f2)

clrs <- yaml::read_yaml("/home/uhlitzf/spectrum_tme/_data/small/signatures/hgsc_v6_colors.yaml") %>% 
  lapply(function(x) map_depth(x, vec_depth(x)-2, unlist))

clrs$patient_id_short <- clrs$patient_id
names(clrs$patient_id_short) <- str_remove_all(names(clrs$patient_id), "SPECTRUM-OV-")
## load meta data
meta_tbl <- read_excel("/home/uhlitzf/spectrum_tme/_data/small/MSK SPECTRUM - Single cell RNA-seq_v6.xlsx", sheet = 2) %>% 
  filter(!(patient_id %in% c("SPECTRUM-OV-100", "SPECTRUM-OV-099")),
         therapy == "pre-Rx") %>% 
  rename(sample = isabl_id) %>% 
  mutate(patient_id_short = str_remove_all(patient_id, "SPECTRUM-OV-"),
         sort_short = str_remove_all(sort_parameters, "singlet, live, "),
         tumor_supersite = str_replace_all(tumor_supersite, "Upper Quadrant", "UQ")) %>% 
  mutate(tumor_supersite = ordered(tumor_supersite, levels = names(clrs$tumor_supersite)))

## load wgs meta data
signature_tbl <- read_tsv("/home/uhlitzf/spectrum_tme/_data/small/mutational_signatures_summary.tsv") %>%
  select(patient_id, consensus_signature) %>% 
  bind_rows(tibble(patient_id = unique(sort(meta_tbl$patient_id[!(meta_tbl$patient_id %in% .$patient_id)])), consensus_signature = "NA")) %>% 
  mutate(consensus_signature = ordered(consensus_signature, levels = names(clrs$consensus_signature))) %>% 
  arrange(consensus_signature) %>% 
  distinct(patient_id, .keep_all = T)

meta_tbl <- left_join(meta_tbl, signature_tbl, by = "patient_id")
## load full seurat objects with expression data
# seu_obj_tc <- read_rds("/work/shah/uhlitzf/data/SPECTRUM/freeze/v5/T.cell_processed_filtered_sub.rds")
# seu_obj_cc <- read_rds(paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v5/Ovarian.cancer.cell_processed_filtered.rds"))
seu_obj_cc <- read_rds("/work/shah/uhlitzf/data/SPECTRUM/freeze/v6/Ovarian.cancer.cell_processed_filtered.rds")
markers_tbl <- read_tsv("/work/shah/uhlitzf/data/SPECTRUM/freeze/v6/Ovarian.cancer.cell_highqc_markers_02.tsv")
plot_data <- cbind(cell_id = colnames(seu_obj_cc), FetchData(seu_obj_cc, c("umapharmony_1", "umapharmony_2", "umappca_1", "umappca_2", "RNA_snn_res.0.2", "sample", "cluster_label", grep("pathway", colnames(seu_obj_cc@meta.data), value = T)))) %>% 
  as_tibble() %>% 
  left_join(meta_tbl, by = "sample") %>% 
  filter(!is.na(consensus_signature))

base_umap <- ggplot(plot_data) +
  coord_fixed() +
  NoAxes() +
  theme(legend.position = c(0, 1),
        legend.justification = c("left", "top"),
        legend.box.just = "left",
        legend.margin = margin(1, 1, 1, 1),
        #panel.border = element_rect(linetype = 1, color = "black", size = 1),
        legend.text = element_text(size = 14, margin = margin(0, 10, 0, 0)),
        legend.spacing.x = unit(0, "npc"), 
        legend.title = element_blank(),
        plot.title = element_text(hjust = 0.5, vjust = 0.5, face = "plain", size = 22))

arrow <- arrow(angle = 20, type = "closed", length = unit(0.1, "npc"))
umap_coord_anno <- ggplot(tibble(group = c("UMAP1", "UMAP2"),
                                 x = c(0, 0), xend = c(1, 0),
                                 y = c(0, 0), yend = c(0, 1),
                                 lx = c(0.5, -0.15), ly = c(-0.15, 0.5),
                                 angle = c(0, 90))) +
  geom_segment(aes(x, y, xend = xend, yend = yend, group = group),
               arrow = arrow, size = 1, lineend = "round") +
  geom_text(aes(lx, ly, label = group, angle = angle), size = 4) +
  theme_void() +
  coord_fixed(xlim = c(-0.3, 1), ylim = c(-0.3, 1))

add_umap_coord <- function(gg_obj) {
  p <- ggdraw() + 
    draw_plot(gg_obj, x = 0, y = 0, width = 1, height = 1) +
    draw_plot(umap_coord_anno, x = -0.015, y = -0.02, width = 0.2, height = 0.2)
  return(p)
}
pt.size <- 0.1
pt.size.mini <- 0.01
pt.alpha <- 0.05
pt.alpha.mini <- 0.02

median_tbl <- plot_data %>%
  group_by(patient_id_short) %>% 
  summarise(umappca_1 = median(umappca_1), 
            umappca_2 = median(umappca_2))

umap_pca_mutsig <- base_umap + 
  geom_point(aes(umappca_1, umappca_2, color = consensus_signature), 
             size = pt.size, alpha = pt.alpha) +
  geom_text(aes(umappca_1, umappca_2, label = patient_id_short), data = median_tbl) + 
  scale_color_manual(values = clrs$consensus_signature) +
  guides(color = guide_legend(override.aes = list(size = 4, alpha = 1),
                              ncol = 1, 
                              label.position = "right")) +
  labs(title = "Mutational signature")

umap_mutsig <- base_umap + 
  geom_point(aes(umapharmony_1, umapharmony_2, color = consensus_signature), 
             size = pt.size, alpha = pt.alpha) +
  scale_color_manual(values = clrs$consensus_signature) +
  guides(color = guide_legend(override.aes = list(size = 4, alpha = 1),
                              ncol = 1, 
                              label.position = "right")) +
  guides(color = F) + 
  labs(title = "")

umap_pca_cluster <- base_umap + 
  geom_point(aes(umappca_1, umappca_2, color = cluster_label), 
             size = pt.size, alpha = pt.alpha) +
  geom_text(aes(umappca_1, umappca_2, label = patient_id_short), data = median_tbl) + 
  scale_color_manual(values = clrs$cluster_label$Ovarian.cancer.cell) +
  guides(color = guide_legend(override.aes = list(size = 4, alpha = 1),
                              ncol = 1, 
                              label.position = "right")) +
  labs(title = "Cluster")

umap_cluster <- base_umap + 
  geom_point(aes(umapharmony_1, umapharmony_2, color = cluster_label), 
             size = pt.size, alpha = pt.alpha) +
  scale_color_manual(values = clrs$cluster_label$Ovarian.cancer.cell) +
  guides(color = guide_legend(override.aes = list(size = 4, alpha = 1),
                              ncol = 1, 
                              label.position = "right")) +
  guides(color = F) + 
  labs(title = "")

cluster_legend <- cowplot::get_legend(umap_pca_cluster)

umap_pca_jak_stat <- base_umap + 
  # geom_point(aes(umappca_1, umappca_2), color = "grey80",
  #            size = pt.size, alpha = pt.alpha, 
  #            data = filter(plot_data, JAK.STAT.pathway <= 0)) +
  geom_point(aes(umappca_1, umappca_2, color = JAK.STAT.pathway), 
             size = pt.size, alpha = pt.alpha, 
             data = mutate(plot_data, JAK.STAT.pathway = ifelse(JAK.STAT.pathway > 4, 4, JAK.STAT.pathway))) +
  geom_text(aes(umappca_1, umappca_2, label = patient_id_short), data = median_tbl) + 
  scale_color_gradientn(colours = viridis(9), breaks = c(0, 2, 4), limits = c(min(plot_data$JAK.STAT.pathway), 4), labels = c(0, 2, "≥4")) +
  labs(title = "JAK-STAT signaling")

umap_jak_stat <- base_umap + 
  # geom_point(aes(umapharmony_1, umapharmony_2), color = "grey80",
  #            size = pt.size, alpha = pt.alpha, 
  #            data = filter(plot_data, JAK.STAT.pathway <= 0)) +
  geom_point(aes(umapharmony_1, umapharmony_2, color = JAK.STAT.pathway), 
             size = pt.size, alpha = pt.alpha, 
             data = mutate(plot_data, JAK.STAT.pathway = ifelse(JAK.STAT.pathway > 4, 4, JAK.STAT.pathway))) +
  scale_color_gradientn(colours = viridis(9), breaks = c(0, 2, 4), limits = c(min(plot_data$JAK.STAT.pathway), 4), labels = c(0, 2, "≥4")) +
  guides(color = F) + 
  labs(title = "")
comp_plot_wrapper <- function(comp_tbl = comp_tbl_sample_sort, cts, y = "nrel", x = "tumor_subsite", fill = "cell_type", switch = NULL, facet = F, exclude_ascites = F, super_type = NULL) {
  if (facet == F) comp_tbl <- filter(comp_tbl, !is.na(tumor_supersite))
  if (exclude_ascites == T) comp_tbl <- filter(comp_tbl, tumor_supersite != "Ascites")
  if (y == "nrel") ylab <- paste0("% cells")
  if (y == "log10n") ylab <- paste0("# cells")
  if (y == "n") ylab <- paste0("# cells")
  p <- ggplot(filter(comp_tbl, sort_short_x == cts),
              aes_string(x, y, fill = fill)) +
    theme(axis.text.x = element_blank(),
          axis.ticks.x = element_blank(),
          axis.title.x = element_blank(),
          axis.line.x = element_blank(),
          axis.title.y = element_text(margin = margin(0, -1, 0, 1, unit = "npc")),
          strip.text.y = element_blank(),
          strip.background.y = element_blank(),
          plot.margin = grid::unit(c(0.04, 0.01, 0.02, 0.01), "npc")) +
    #scale_x_discrete(expand = c(0, 0)) +
    labs(x = "",
         y = ylab) +
    guides(fill = FALSE)
  if (fill != "cluster_label") p <- p + scale_fill_manual(values = clrs[[fill]])
  if (fill == "cluster_label") p <- p + scale_fill_manual(values = clrs[[fill]][[super_type]])
  if (facet == "tumor_supersite") p <- p +
    facet_grid(sort_short_x~tumor_supersite,
               space = "free", scales = "free", switch = switch)
  if (facet == "consensus_signature") p <- p +
    facet_grid(sort_short_x~consensus_signature,
               space = "free", scales = "free", switch = switch)
  if (facet == "patient_id_short") p <- p +
    facet_grid(sort_short_x~patient_id_short,
               space = "free", scales = "free", switch = switch)
  if (facet == F) p <- p +
    facet_grid(sort_short_x~.,
               space = "free", scales = "free", switch = switch)
  if (y == "nrel") p <- p +
    geom_bar(stat = "identity", width = 1) +
    scale_y_continuous(expand = c(0, 0),
                       breaks = c(0, 50, 100),
                       labels = c("0", "50", "100")) +
    theme(strip.text.x = element_blank(),
          strip.background.x = element_blank())
  if (y == "log10n") p <- p +
    geom_bar(stat = "identity", position = position_dodge2(preserve = "single"), width = 1) +
    scale_y_continuous(expand = c(0, 0),
                       breaks = c(1, 2, 3),
                       limits = c(0, 4),
                       labels = function(x) parse(text = paste("10^", x))) +
    expand_limits(y = c(0, 4)) +
    theme(panel.grid.major.y = element_line(linetype = 1, color = "grey90", size = 0.5))
  if (y == "n") p <- p +
    geom_bar(stat = "identity", position = position_stack(), width = 1) +
    scale_y_continuous(expand = c(0, 0),
                       breaks = c(0, 2500, 5000, 7500),
                       limits = c(0, 7500),
                       labels = c("", "2500", "5000", "7500")) +
    expand_limits(y = c(0, 7500)) +
    theme(panel.grid.major.y = element_line(linetype = 1, color = "grey90", size = 0.5),
          plot.title = element_text(face = "plain", size = 14,
                                    hjust = 0.5, vjust = 0.5)) +
    labs(title = cts)
  return(p)
}

common_label_layers <- list(
  theme(axis.title.y = element_blank(),
        strip.text = element_blank(),
        strip.background = element_blank(),
        axis.text.x = element_blank(),
        axis.title.x = element_blank(),
        axis.ticks = element_blank(),
        axis.line.x = element_blank(),
        plot.margin = ggplot2::margin(0.02, 0.01, 0.02, 0.01, "npc")),
  scale_y_discrete(expand = c(0, 0)),
  guides(fill = F, color = F)
)

comp_label <- function(x = "tumor_subsite", y = "label_supersite", fill = "tumor_supersite", facet = F, cts = c("CD45-", "CD45+")) {
  comp_tbl <- comp_tbl_sample_sort
  if (facet == F) comp_tbl <- filter(comp_tbl, !is.na(tumor_supersite))
  p <- ggplot(distinct(filter(comp_tbl, therapy == "pre-Rx", sort_short_x %in% cts), tumor_subsite, tumor_supersite, consensus_signature, label_supersite, label_mutsig, patient_id_short, sample_id_Myeloid.cell, sample_id_T.cell, rank_T.cell, sample_id_Fibroblast, sample_id_Ovarian.cancer.cell, rank_Ovarian.cancer.cell, sample_id_B.cell, sample_id_Plasma.cell),
              aes_string(x, y, fill = fill)) +
    geom_tile() +
    common_label_layers
  if (fill != "cluster_label") p <- p + scale_fill_manual(values = clrs[[fill]])
  if (fill == "cluster_label") p <- p + scale_fill_manual(values = clrs[[fill]][[super_type]])
  if (facet == "patient_id_short") p <- p + facet_grid(label_supersite~patient_id_short,
                                      space = "free", scales = "free")
  if (facet == "tumor_supersite") p <- p + facet_grid(label_supersite~tumor_supersite,
                                      space = "free", scales = "free")
  if (facet == F) p <- p + facet_grid(label_supersite~.,
                                      space = "free", scales = "free")
  return(p)
}

comp_label_boxplot <- function(comp_tbl = comp_tbl_sample_sort, x = "tumor_supersite", y = "rank_T.cell", facet = T, cts = c("CD45-", "CD45+"), clr = "tumor_supersite", ct = names(clrs$cell_type), pvals = F, exclude_ascites = F) {
  if (facet == F) comp_tbl <- filter(comp_tbl, !is.na(tumor_supersite))
  if (exclude_ascites == T) comp_tbl <- filter(comp_tbl, tumor_supersite != "Ascites")
  # rank_tests_sub <- filter(rank_tests, sort_short_x %in% cts, pval < 0.05) %>%
  #   mutate(pstar = "*")
  required_columns <- c("tumor_subsite", "tumor_supersite", "consensus_signature", "BRCA1_status", "BRCA2_status", "BRCA_status", "label_supersite", "patient_id_short", "sample_id_Myeloid.cell", "sample_id_T.cell", "rank_T.cell", "sample_id_Fibroblast", "sample_id_Ovarian.cancer.cell", "rank_Ovarian.cancer.cell", "sample_id_B.cell", "sample_id_Plasma.cell", "rank_IR", "sample_id_IR", "rank_MC", "sample_id_MC", "sample_cd45p", "sample_rank")
  missing_columns <- required_columns[!(required_columns %in% colnames(comp_tbl))]
  comp_tbl <- bind_cols(comp_tbl, matrix(rep("distinct_helper", n = length(missing_columns)), ncol = length(missing_columns)) %>% set_colnames(missing_columns) %>% as_tibble())
  p <- ggplot(
    distinct(filter(comp_tbl, therapy == "pre-Rx", sort_short_x %in% cts, cell_type %in% ct), tumor_subsite, tumor_supersite, consensus_signature, BRCA1_status, BRCA2_status, BRCA_status, label_supersite, patient_id_short, sample_id_Myeloid.cell, sample_id_T.cell, rank_T.cell, sample_id_Fibroblast, sample_id_Ovarian.cancer.cell, rank_Ovarian.cancer.cell, sample_id_B.cell, sample_id_Plasma.cell, rank_IR, sample_id_IR, rank_MC, sample_id_MC, sample_cd45p, sample_rank)
     # distinct(filter(comp_tbl, therapy == "pre-Rx", sort_short_x %in% cts, cell_type %in% ct), -n, -log10n, -nrel)
    ) +
    geom_tile(aes_string(x, y, fill = clr),
              alpha = 0.8) +
    geom_boxplot(aes_string(x, y, fill = clr),
                 width = 0.35, size = 2.5, color = "white", outlier.shape = NA) +
    geom_boxplot(aes_string(x, y, color = clr),
                 width = 0.35, size = 1, fill = "white", outlier.shape = NA) +
    geom_hline(yintercept = 0, linetype = 2, size = 0.5, color = "black") +
    scale_fill_manual(values = clrs[[clr]]) +
    scale_color_manual(values = clrs[[clr]]) +
    coord_flip(clip = "off", ylim = c(-1.01, 1.01)) +
    theme(axis.title.y = element_blank(),
          strip.text = element_blank(),
          strip.background = element_blank(),
          plot.margin = ggplot2::margin(0.04, 0.05, 0.02, 0.01, "npc")) +
    scale_x_discrete(expand = c(0, 0)) +
    scale_y_continuous(expand = c(0, 0), breaks = c(-1, 0, 1)) +
    guides(fill = F, color = F) +
    labs(y = "Scaled rank")
  if (facet == T) p <- p + facet_grid(label_supersite~patient_id_short,
                                      space = "free", scales = "free")
  if (facet == F) p <- p + facet_grid(label_supersite~.,
                                      space = "free", scales = "free")
  if (pvals == T) p <- p +
    geom_text(aes(x = comparison, y = 1.05, label = pstar), data = rank_tests_sub, hjust = 0.5, vjust = 1, size = 5, angle = 90) #+
    # geom_segment(aes(x = "Adnexa", xend = comparison, y = 1.07, yend = 1.07), data = rank_tests_sub)
  return(p)
}


# ggdraw() + draw_plot(comp_label_boxplot(x = "tumor_supersite", y = "rank_T.cell", facet = F, cts = "CD45+", ct = "T cell", pvals = T), width = 0.9)

hist_plot_wrapper <- function(cts, x = "nrel") {
  if (x == "nrel") {
    xlab <- paste0("% cells")
    bw <- 5
  }
  if (x == "log10n") {
    xlab <- paste0("# cells")
    bw <- 0.2
  }
  p <- ggplot(filter(comp_tbl_sample_sort, sort_short_x == cts, !is.na(cell_type))) +
    ggridges::geom_density_ridges(
      aes_string(x, "cell_type", fill = "cell_type"), color = "black",
      stat = "binline", binwidth = bw, scale = 3) +
    facet_grid(label_supersite~label_therapy,
               space = "free", scales = "free") +
    scale_fill_manual(values = clrs$cell_type) +
    theme(axis.text.y = element_blank(),
          axis.ticks.y = element_blank(),
          axis.title.y = element_blank(),
          axis.line.y = element_blank(),
          strip.text = element_blank(),
          strip.background = element_blank(),
          plot.margin = ggplot2::margin(0.02, 0, 0.02, 0, "npc")) +
    scale_x_continuous(expand = c(0, 0)) +
    scale_y_discrete(expand = c(0, 0)) +
    guides(fill = FALSE) +
    labs(x = xlab)
  if (x == "log10n") p <- p + expand_limits(x = c(0, 4)) +
    scale_x_continuous(expand = c(0, 0),
                       labels = function(x) parse(text = paste("10^", x)))
  return(p)
}
cluster_comp <- plot_data %>% 
  mutate(sort_short_x = str_replace_all(sort_short, "U", "CD45-")) %>% 
  group_by(cluster_label, consensus_signature, sample, sort_short_x, tumor_supersite, therapy) %>% 
  tally %>% 
  group_by(consensus_signature, sample, sort_short_x, tumor_supersite, therapy) %>% 
  mutate(nrel = n/sum(n)*100) %>% 
  ungroup %>% 
  arrange(cluster_label != "Secretory.cell.6", -nrel) %>% 
  mutate(sample = ordered(sample, levels = unique(sample))) %>% 
  mutate(cell_type = "Ovarian.cancer.cell")

cluster_comp_rank <- cluster_comp %>% 
  distinct(sample, .keep_all = T) %>% 
  filter(sort_short_x == "CD45-") %>% 
  mutate(sample = ordered(sample, levels = unique(as.character(sample)))) %>% 
  mutate(sample_rank = scales::rescale(as.numeric(sample), c(-1, 1)))

p1 <- comp_plot_wrapper(cluster_comp, "CD45-", "n", x = "sample", fill = "cluster_label", super_type = "Ovarian.cancer.cell")
p2 <- comp_plot_wrapper(cluster_comp, "CD45-", "nrel", x = "sample", fill = "cluster_label", super_type = "Ovarian.cancer.cell")
p3 <- comp_label_boxplot(cluster_comp_rank, x = "consensus_signature", y = "sample_rank", clr = "consensus_signature", cts = "CD45-")


pcomp_grid_p_full <- plot_grid(p1, p2, p3,
                               ncol = 1, align = "v",
                               rel_heights = c(0.25, 0.25, 0.4))
# cluster_comp <- plot_data %>% 
#   group_by(cluster_label, consensus_signature) %>% 
#   tally %>% 
#   group_by(consensus_signature) %>% 
#   mutate(total_cs = sum(n)) %>% 
#   group_by(cluster_label) %>% 
#   # filter(total_cs > 10000) %>% 
#   mutate(total_cl = sum(n),
#          nnorm = n/total_cs,
#          nrel = n/sum(n),
#          nrelnorm = nnorm/sum(nnorm)) %>% 
#   gather(key, value, -cluster_label, -consensus_signature, -total_cs, -total_cl)
# 
# cluster_comp %>% 
#   filter(key %in% c("nrelnorm"),
#          total_cl > 10000) %>% 
#   ggplot() +
#   geom_bar(aes(cluster_label, value, fill = consensus_signature), 
#            stat = "identity", position = "dodge") +
#   facet_wrap(~key, scales = "free_x") + 
#   scale_fill_manual(values = clrs$consensus_signature) +
#   theme(axis.ticks.y = element_blank(),
#         axis.title.x = element_text(),
#         axis.text.y = element_text(hjust = 1)) +
#   coord_flip()
cancer_grid_pca <- ggdraw() +
  draw_plot(add_umap_coord(umap_pca_mutsig), 
            x = 0, y = 0, width = 0.25, height = 1) +
  draw_plot(add_umap_coord(umap_pca_cluster + guides(color = F)), 
            x = 0.25, y = 0, width = 0.25, height = 1) +
  draw_grob(cluster_legend, x = 0.5, y = -0.15, height = 1) +
  draw_plot(add_umap_coord(umap_pca_jak_stat), x = 0.7, y = 0, width = 0.25, height = 1)

cancer_grid <- ggdraw() +
  draw_plot(add_umap_coord(umap_mutsig), 
            x = 0, y = 0, width = 0.25, height = 1) +
  draw_plot(add_umap_coord(umap_cluster + guides(color = F)), 
            x = 0.25, y = 0, width = 0.25, height = 1) +
  draw_plot(pcomp_grid_p_full, 
            x = 0.5, y = 0, width = 0.18, height = 1) +
  draw_plot(add_umap_coord(umap_jak_stat), x = 0.7, y = 0, width = 0.25, height = 1)

cancer_grid_pca

cancer_grid

# ggsave("_fig/003_cancer_cell/003_umap_grid.pdf", cancer_grid, width = 20, height = 5)
ggsave("_fig/003_cancer_cell/003_umap_grid_pca.png", cancer_grid_pca, width = 20, height = 5)
ggsave("_fig/003_cancer_cell/003_umap_grid.png", cancer_grid, width = 20, height = 5)
# ## downsample objects for exploration
# # set.seed(123)
# # seu_obj_tc_sub <- subset(seu_obj_tc, cells = sample(colnames(seu_obj_tc), 10000))
# # seu_obj_cc_sub <- subset(seu_obj_cc, cells = sample(colnames(seu_obj_cc), 10000))
# seu_obj_cc_sub <- seu_obj_cc
# seu_obj_tc_sub <- seu_obj_tc
# 
# tc_data <- cbind(cell_id = colnames(seu_obj_tc_sub), FetchData(seu_obj_tc_sub, c("sample", "cluster_label", names(signature_modules), "CD8.ISG"))) %>% 
#   filter(str_detect(cluster_label, "CD8")) %>% 
#   as_tibble %>% 
#   rename(isabl_id = sample) %>% 
#   gather(module, tc_score, -cell_id, -isabl_id, -cluster_label) %>% 
#   left_join(meta_tbl, by = "isabl_id")
# 
# patient_tc_ranks <- tc_data %>% 
#   filter(sort_short == "CD45+", therapy == "pre-Rx") %>% 
#   group_by(sample_id, module) %>% 
#   summarise(median_tc_score = median(tc_score)) %>% 
#   ungroup

1 Differential expression in cancer cells

Are there certain pathways activated in cancer cells…

# progeny_tbl <- progeny(as.matrix(seu_obj_cc_sub@assays$RNA@data[VariableFeatures(seu_obj_cc_sub),])) %>% 
#   as_tibble(rownames = "cell_id") %>% 
#   gather(pathway, score, -cell_id)
# 
# plot_data <- cbind(cell_id = colnames(seu_obj_cc_sub), FetchData(seu_obj_cc_sub, c("umapharmony_1", "umapharmony_2", "sample", "RNA_snn_res.0.2", "nFeature_RNA", "percent.mt"))) %>% 
#   as_tibble() %>% 
#   rename(isabl_id = sample) %>% 
#   left_join(meta_tbl, by = "isabl_id") %>% 
#   filter(sort_short == "CD45-", therapy == "pre-Rx") %>% 
#   left_join(select(signature_tbl, patient_id, consensus_signature), by = "patient_id") %>% 
#   left_join(progeny_tbl, by = "cell_id")
# 
# patient_pw_ranks <- plot_data %>% 
#   group_by(pathway, patient_id) %>% 
#   summarise(median_score = median(score)) %>% 
#   ungroup() %>% 
#   arrange(pathway, median_score) %>% 
#   # group_by(pathway) %>% 
#   mutate(pw_rank = as.factor(row_number(median_score))) %>% 
#   ungroup()
# 
# patient_pw_ranks_site <- plot_data %>% 
#   group_by(pathway, sample_id) %>% 
#   summarise(median_score_site = median(score)) %>% 
#   ungroup() %>% 
#   arrange(pathway, median_score_site) %>% 
#   # group_by(pathway) %>% 
#   mutate(pw_rank_site = as.factor(row_number(median_score_site))) %>% 
#   ungroup()
# 
# plot_data_ranked <- plot_data %>% 
#   left_join(patient_pw_ranks, by = c("pathway", "patient_id")) %>% 
#   left_join(patient_pw_ranks_site, by = c("pathway", "sample_id"))

1.1 pathway scores across cancer cell clusters

# p1 <- ggplot() + 
#   geom_point(aes(umapharmony_1, umapharmony_2, color = RNA_snn_res.0.2), size = 0.01, 
#              data = distinct(plot_data, umapharmony_1, umapharmony_2, RNA_snn_res.0.2)) +
#   theme_void() +
#   guides(color = guide_legend(override.aes = list(size = 2))) +
#   coord_fixed()
# 
# p2 <- ggplot() +
#   geom_bar(aes(RNA_snn_res.0.2, n, fill = RNA_snn_res.0.2), stat = "identity",
#            data = distinct(plot_data, umapharmony_1, umapharmony_2, RNA_snn_res.0.2) %>% 
#              group_by(RNA_snn_res.0.2) %>% tally) +
#   theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.5)) +
#   scale_y_log10() + 
#   guides(fill = F)
# 
# p3 <- ggplot() + 
#   geom_point(aes(umapharmony_1, umapharmony_2), size = 0.01, color = "grey80",
#              data = distinct(plot_data, umapharmony_1, umapharmony_2)) +
#   geom_point(aes(umapharmony_1, umapharmony_2, color = RNA_snn_res.0.2), size = 0.01, alpha = 0.1,
#              data = distinct(plot_data, umapharmony_1, umapharmony_2, RNA_snn_res.0.2)) +
#   facet_wrap(~RNA_snn_res.0.2, ncol = 8) + 
#   theme_void() +
#   guides(color = guide_legend(override.aes = list(size = 2))) +
#   coord_fixed()
# 
# plot_grid(p1, p2, nrow = 1)
# p3
# ggplot() + 
#   geom_point(aes(umapharmony_1, umapharmony_2), size = 0.1, color = "grey80",
#              data = filter(plot_data, score <= 0)) +
#   geom_point(aes(umapharmony_1, umapharmony_2, color = score), size = 0.1, 
#              data = filter(plot_data, score > 0)) +
#   facet_wrap(~pathway) + 
#   scale_color_gradientn(colors = viridis(9)) + 
#   theme_void()
# 
# ggplot() + 
#   geom_boxplot(aes(RNA_snn_res.0.2, score, color = RNA_snn_res.0.2), 
#                data = plot_data) +
#   facet_wrap(~pathway, scales = "free") + 
#   # scale_color_gradientn(colors = viridis(9)) + 
#   theme_void() +
#   theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.5))
# 
# select(markers_tbl, cluster, gene) %>% 
#   group_by(cluster) %>% 
#   mutate(rank = row_number(cluster)) %>% 
#   slice(1:20) %>% 
#   spread(cluster, gene) %>% 
#   knitr::kable()

1.2 pathway scores across tumor sites

# ggplot() + 
#   geom_boxplot(aes(tumor_supersite, score, color = tumor_supersite),
#                data = plot_data_ranked) +
#   facet_wrap(~pathway, scales = "free_y") + 
#   scale_color_manual(values = clrs$tumor_supersite) + 
#   theme_void()
# 
# ggplot() + 
#   geom_boxplot(aes(pw_rank_site, score, color = tumor_supersite, group = isabl_id), 
#                data = plot_data_ranked) +
#   facet_wrap(~pathway, scales = "free") + 
#   scale_color_manual(values = clrs$tumor_supersite) + 
#   theme_void()

1.3 pathway scores across mutational signatures

# ggplot() + 
#   geom_boxplot(aes(consensus_signature, score, color = consensus_signature), 
#                data = plot_data) +
#   facet_wrap(~pathway, scales = "free_y") + 
#   scale_color_manual(values = clrs$consensus_signature) + 
#   theme_void()
# 
# ggplot(plot_data_ranked) + 
#   geom_boxplot(aes(pw_rank, score, color = consensus_signature, group = patient_id)) +
#   facet_wrap(~pathway, scales = "free") + 
#   scale_color_manual(values = clrs$consensus_signature) + 
#   scale_x_discrete(breaks = arrange(distinct(plot_data_ranked, pw_rank, patient_id_short), pw_rank)$pw_rank,
#                      labels = arrange(distinct(plot_data_ranked, pw_rank, patient_id_short), pw_rank)$patient_id_short) +
#   theme_void() +
#   theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.5))

1.4 correlation of pathway scores and T cell scores

# tc_cc_cor_tbl <- patient_pw_ranks_site %>% 
#   left_join(patient_tc_ranks, by = "sample_id") %>% 
#   left_join(distinct(meta_tbl, patient_id, sample_id, tumor_supersite), 
#             by = "sample_id") %>% 
#   left_join(select(signature_tbl, patient_id, consensus_signature), 
#             by = "patient_id") %>% 
#   na.omit()
# 
# ggplot(tc_cc_cor_tbl) +
#   geom_point(aes(median_tc_score, median_score_site, color = tumor_supersite),
#              size = 0.1) +
#   facet_grid(module~pathway) +
#   scale_color_manual(values = clrs$tumor_supersite) +
#   stat_cor(aes(median_tc_score, median_score_site), method = "spearman",
#            label.sep = "\n") +
#   geom_smooth(aes(median_tc_score, median_score_site), 
#               method = "lm", color = "black")
# 
# ggplot(tc_cc_cor_tbl) +
#   geom_point(aes(median_tc_score, median_score_site, color = consensus_signature),
#              size = 0.1) +
#   facet_grid(module~pathway) +
#   scale_color_manual(values = clrs$consensus_signature) +
#   stat_cor(aes(median_tc_score, median_score_site), method = "spearman",
#            label.sep = "\n") +
#   geom_smooth(aes(median_tc_score, median_score_site), 
#               method = "lm", color = "black")

2 JAK-STAT signaling

# p1 <- ggplot() + 
#   geom_boxplot(aes(consensus_signature, score, color = consensus_signature), 
#                data = filter(plot_data, pathway == "JAK-STAT")) +
#   facet_wrap(~pathway, scales = "free_y") + 
#   scale_color_manual(values = clrs$consensus_signature) +
#   theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.5)) +
#   guides(color = F) +
#   labs(x = "", y = "PROGENy score")
# 
# plot_data_ranked_jak_stat <- filter(plot_data_ranked, pathway == "JAK-STAT")
# 
# p2 <- ggplot(plot_data_ranked_jak_stat) + 
#   geom_boxplot(aes(pw_rank, score, color = consensus_signature, group = patient_id)) +
#   facet_wrap(~pathway, scales = "free") + 
#   scale_color_manual(values = clrs$consensus_signature) + 
#   scale_x_discrete(breaks = arrange(distinct(plot_data_ranked_jak_stat, pw_rank, patient_id_short), pw_rank)$pw_rank,
#                      labels = arrange(distinct(plot_data_ranked_jak_stat, pw_rank, patient_id_short), pw_rank)$patient_id_short) +
#   theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.5)) +
#   guides(color = F) +
#   labs(x = "", y = "PROGENy score", color = "Signature")
# 
# p3 <- ggplot(filter(tc_cc_cor_tbl, pathway == "JAK-STAT", module == "CD8.ISG")) +
#   geom_point(aes(median_tc_score, median_score_site, color = consensus_signature),
#              size = 1) +
#   facet_wrap(~pathway) +
#   scale_color_manual(values = clrs$consensus_signature) +
#   stat_cor(aes(median_tc_score, median_score_site), method = "spearman",
#            label.sep = ", ") +
#   geom_smooth(aes(median_tc_score, median_score_site), 
#               method = "lm", color = "black") +
#   labs(x = "ISG module score\nin CD8 T cells", 
#        y = "JAK-STAT PROGENy\nscore in cancer cells", 
#        color = "Signature") +
#   theme(strip.text = element_blank())
# 
# plot_grid(plot_grid(p1, p2, align = "h", rel_widths = c(0.3, 0.7)), 
#           p3, rel_widths = c(0.6, 0.4))

3 session info

devtools::session_info()
## ─ Session info ───────────────────────────────────────────────────────────────
##  setting  value                       
##  version  R version 3.6.2 (2019-12-12)
##  os       Debian GNU/Linux 10 (buster)
##  system   x86_64, linux-gnu           
##  ui       X11                         
##  language (EN)                        
##  collate  en_US.UTF-8                 
##  ctype    en_US.UTF-8                 
##  tz       Etc/UTC                     
##  date     2020-12-11                  
## 
## ─ Packages ───────────────────────────────────────────────────────────────────
##  ! package      * version    date       lib
##    abind          1.4-5      2016-07-21 [2]
##    ape            5.3        2019-03-17 [2]
##    assertthat     0.2.1      2019-03-21 [2]
##    backports      1.1.10     2020-09-15 [1]
##    bibtex         0.4.2.2    2020-01-02 [2]
##    Biobase        2.46.0     2019-10-29 [2]
##    BiocGenerics   0.32.0     2019-10-29 [2]
##    bitops         1.0-6      2013-08-17 [2]
##    broom          0.7.2      2020-10-20 [1]
##    callr          3.4.2      2020-02-12 [1]
##    car            3.0-8      2020-05-21 [1]
##    carData        3.0-4      2020-05-22 [1]
##    caTools        1.17.1.4   2020-01-13 [2]
##    cellranger     1.1.0      2016-07-27 [2]
##    cli            2.0.2      2020-02-28 [1]
##    cluster        2.1.0      2019-06-19 [3]
##    codetools      0.2-16     2018-12-24 [3]
##    colorblindr  * 0.1.0      2020-01-13 [2]
##    colorspace   * 1.4-2      2019-12-29 [2]
##    cowplot      * 1.0.0      2019-07-11 [2]
##    crayon         1.3.4      2017-09-16 [1]
##    curl           4.3        2019-12-02 [2]
##    data.table     1.12.8     2019-12-09 [2]
##    DBI            1.1.0      2019-12-15 [2]
##    dbplyr         2.0.0      2020-11-03 [1]
##    desc           1.2.0      2018-05-01 [2]
##    devtools       2.2.1      2019-09-24 [2]
##    digest         0.6.25     2020-02-23 [1]
##    dplyr        * 1.0.2      2020-08-18 [1]
##    ellipsis       0.3.1      2020-05-15 [1]
##    evaluate       0.14       2019-05-28 [2]
##    fansi          0.4.1      2020-01-08 [2]
##    farver         2.0.3      2020-01-16 [1]
##    fitdistrplus   1.0-14     2019-01-23 [2]
##    forcats      * 0.5.0      2020-03-01 [1]
##    foreign        0.8-74     2019-12-26 [3]
##    fs             1.5.0      2020-07-31 [1]
##    future         1.15.1     2019-11-25 [2]
##    future.apply   1.4.0      2020-01-07 [2]
##    gbRd           0.4-11     2012-10-01 [2]
##    gdata          2.18.0     2017-06-06 [2]
##    generics       0.0.2      2018-11-29 [2]
##    ggplot2      * 3.3.2      2020-06-19 [1]
##    ggpubr       * 0.4.0      2020-06-27 [1]
##    ggrepel        0.8.1      2019-05-07 [2]
##    ggridges       0.5.2      2020-01-12 [2]
##    ggsignif       0.6.0      2019-08-08 [1]
##    globals        0.12.5     2019-12-07 [2]
##    glue           1.3.2      2020-03-12 [1]
##    gplots         3.0.1.2    2020-01-11 [2]
##    gridExtra      2.3        2017-09-09 [2]
##    gtable         0.3.0      2019-03-25 [2]
##    gtools         3.8.1      2018-06-26 [2]
##    haven          2.3.1      2020-06-01 [1]
##    hms            0.5.3      2020-01-08 [1]
##    htmltools      0.4.0      2019-10-04 [2]
##    htmlwidgets    1.5.1      2019-10-08 [2]
##    httr           1.4.2      2020-07-20 [1]
##    ica            1.0-2      2018-05-24 [2]
##    igraph         1.2.5      2020-03-19 [1]
##    irlba          2.3.3      2019-02-05 [2]
##    jsonlite       1.7.1      2020-09-07 [1]
##    KernSmooth     2.23-16    2019-10-15 [3]
##    knitr          1.26       2019-11-12 [2]
##    labeling       0.3        2014-08-23 [2]
##    lattice        0.20-38    2018-11-04 [3]
##    lazyeval       0.2.2      2019-03-15 [2]
##    leiden         0.3.1      2019-07-23 [2]
##    lifecycle      0.2.0      2020-03-06 [1]
##    listenv        0.8.0      2019-12-05 [2]
##    lmtest         0.9-37     2019-04-30 [2]
##    lsei           1.2-0      2017-10-23 [2]
##    lubridate      1.7.9.2    2020-11-13 [1]
##  P magick       * 2.4.0      2020-06-23 [?]
##    magrittr     * 2.0.1      2020-11-17 [1]
##    MASS           7.3-51.5   2019-12-20 [3]
##    Matrix         1.2-18     2019-11-27 [3]
##    memoise        1.1.0      2017-04-21 [2]
##    metap          1.2        2019-12-08 [2]
##    mnormt         1.5-5      2016-10-15 [2]
##    modelr         0.1.8      2020-05-19 [1]
##    multcomp       1.4-12     2020-01-10 [2]
##    multtest       2.42.0     2019-10-29 [2]
##    munsell        0.5.0      2018-06-12 [2]
##    mutoss         0.1-12     2017-12-04 [2]
##    mvtnorm        1.0-12     2020-01-09 [2]
##    nlme           3.1-143    2019-12-10 [3]
##    npsurv         0.4-0      2017-10-14 [2]
##    numDeriv       2016.8-1.1 2019-06-06 [2]
##    openxlsx       4.1.5      2020-05-06 [1]
##    pbapply        1.4-2      2019-08-31 [2]
##    pillar         1.4.6      2020-07-10 [1]
##    pkgbuild       1.0.6      2019-10-09 [2]
##    pkgconfig      2.0.3      2019-09-22 [1]
##    pkgload        1.0.2      2018-10-29 [2]
##    plotly         4.9.1      2019-11-07 [2]
##    plotrix        3.7-7      2019-12-05 [2]
##    plyr           1.8.5      2019-12-10 [2]
##    png            0.1-7      2013-12-03 [2]
##    prettyunits    1.1.1      2020-01-24 [1]
##    processx       3.4.2      2020-02-09 [1]
##    ps             1.3.2      2020-02-13 [1]
##    purrr        * 0.3.4      2020-04-17 [1]
##    R.methodsS3    1.7.1      2016-02-16 [2]
##    R.oo           1.23.0     2019-11-03 [2]
##    R.utils        2.9.2      2019-12-08 [2]
##    R6             2.4.1      2019-11-12 [1]
##    RANN           2.6.1      2019-01-08 [2]
##    rappdirs       0.3.1      2016-03-28 [2]
##    RColorBrewer   1.1-2      2014-12-07 [2]
##    Rcpp           1.0.4      2020-03-17 [1]
##    RcppAnnoy      0.0.16     2020-03-08 [1]
##    RcppParallel   4.4.4      2019-09-27 [2]
##    Rdpack         0.11-1     2019-12-14 [2]
##    readr        * 1.4.0      2020-10-05 [1]
##    readxl       * 1.3.1      2019-03-13 [2]
##    remotes        2.1.0      2019-06-24 [2]
##    reprex         0.3.0      2019-05-16 [2]
##    reshape2       1.4.3      2017-12-11 [2]
##    reticulate     1.14       2019-12-17 [2]
##    rio            0.5.16     2018-11-26 [1]
##    rlang          0.4.8      2020-10-08 [1]
##    rmarkdown      2.0        2019-12-12 [2]
##    ROCR           1.0-7      2015-03-26 [2]
##    rprojroot      1.3-2      2018-01-03 [2]
##    rstatix        0.6.0      2020-06-18 [1]
##    rstudioapi     0.11       2020-02-07 [1]
##    rsvd           1.0.3      2020-02-17 [1]
##    Rtsne          0.15       2018-11-10 [2]
##    rvest          0.3.6      2020-07-25 [1]
##    sandwich       2.5-1      2019-04-06 [2]
##    scales         1.1.0      2019-11-18 [2]
##    sctransform    0.2.1      2019-12-17 [2]
##    SDMTools       1.1-221.2  2019-11-30 [2]
##    sessioninfo    1.1.1      2018-11-05 [2]
##    Seurat       * 3.1.2      2019-12-12 [2]
##    sn             1.5-4      2019-05-14 [2]
##    stringi        1.5.3      2020-09-09 [1]
##    stringr      * 1.4.0      2019-02-10 [1]
##    survival       3.1-8      2019-12-03 [3]
##    testthat       2.3.2      2020-03-02 [1]
##    TFisher        0.2.0      2018-03-21 [2]
##    TH.data        1.0-10     2019-01-21 [2]
##    tibble       * 3.0.4      2020-10-12 [1]
##    tidyr        * 1.1.2      2020-08-27 [1]
##    tidyselect     1.1.0      2020-05-11 [1]
##    tidyverse    * 1.3.0      2019-11-21 [2]
##    tsne           0.1-3      2016-07-15 [2]
##    usethis        1.5.1      2019-07-04 [2]
##    uwot           0.1.5      2019-12-04 [2]
##    vctrs          0.3.5      2020-11-17 [1]
##    viridis      * 0.5.1      2018-03-29 [2]
##    viridisLite  * 0.3.0      2018-02-01 [2]
##    withr          2.3.0      2020-09-22 [1]
##    xfun           0.12       2020-01-13 [2]
##    xml2           1.3.2      2020-04-23 [1]
##    yaml           2.2.1      2020-02-01 [1]
##    zip            2.0.4      2019-09-01 [1]
##    zoo            1.8-7      2020-01-10 [2]
##  source                                 
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.3)                         
##  CRAN (R 3.6.2)                         
##  Bioconductor                           
##  Bioconductor                           
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  Github (clauswilke/colorblindr@1ac3d4d)
##  R-Forge (R 3.6.2)                      
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.3)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.3)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.3)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.3)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.3)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  Bioconductor                           
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.3)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.3)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.3)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.3)                         
##  CRAN (R 3.6.3)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.3)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
## 
## [1] /home/uhlitzf/R/lib
## [2] /usr/local/lib/R/site-library
## [3] /usr/local/lib/R/library
## 
##  P ── Loaded and on-disk path mismatch.